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. i-O.' INTRODUCTION 

I’uture NASA geo dynamic and oceanographic programs depend on extremely 
accurate knowledge of the geoid. Much of this knowledge will come from highly 
accurate altimeter measurements from future satellite missions ( i*e. , SEASAT) . 

The full use of this data necessitates realistic models of the global 
ocean tides: for both the reduction of altimeter data to mean sea level and , 
for direct orbital perturbations. Moreover, global ocean tide models are of 

V 

interest In the study of energy dissipation in the earth-moon system and, for 
oceanographic studies, the global ocean tidal currents are of interest - 

The dynamics of the global ocean tides are described by the nonlinear 
Laplace Tidal Partial Differential Equations , The equations are too complex 
to solve analytically^ and are amenable only to numerical solution by computer. 

In a previous work (Estes [9]) which we will refer to as report [A], a computer 
software system was developed to solve for theoretical global tides described by 
the Laplace Tidal Equations (LTE) including static earth tides, bottom friction and 
turbulent viscosity. Following Zahel [25], a time stepping integration 
technique incorporating impermeable coastal boundary conditions and a forward 
time difference/ central spatial difference numerical scheme was used. The 
present development is an extension of that work to include ocean self-gravitation 
and crustal loading effects into the Laplace Tidal Equations, As discussed 
by Hendershott jl5] , the inclusion of these effects forces the LTE into an 
integro-differential system which we denote by IDLTE (^ntegro - ^fferential 
' Laplace ^idal Equation) , For purposes of completene^, much of the material of . 
report [A] will be repeated, particularly in the development of the tidal 
equations and in the description of the software system design. 


Report [A] presented numerical solutions of LTE for the M^ and tidal 
constituents* Ne have since obtained solutions for 0^ and For 

convenience of reference, we include a complete set of these theoretical 
solutions, which show general agreement with empirical global solutions with 
respect to the positions of amphidromic systems and tidal amplitudes. An M^ 
solution of IDLTE (i.e., including ocean self-gravitation and crustal loading) 
has been generated and displays a marked improvement in the predicted phase 


relations with observations (generally within 30°) when compared to the LTE 
solution, IDLTE solutions for additional constituents will be generated 


and presented in a future publication- 


The iinprpved comparison with observations obtained from the theoretical 
solution with the inclusion of self “gravitation and crustal loading effects 
are encouraging. However, the observed values of ocean tides contain meteorological 
effects and effects from terms beyond the second order term of the astronomical 
tidal potential which are neglected in the theoretical models. Moreover,- 
internal waves, tide currents and long period ocean tides are not considered 
in the present study. In addition, the parameters representing the earth 
model (Love numbers and load Love numbers) are averaged quantities so that 
local responses of the earth due to inhomogeneities are not included. These 
effects must be correctly modeled before agreement to the 10 centimeter level 
can be achieved. 
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2, 0 LAPI^^CE !TIDAL EQUATIONS 

The dynasiic equations which describe the ocean tides ^ first derived by 
Laplace (1775) , are • 



2wsin6v — ' — — 1 + f-i 

&t ^ Rcos(> ax J 


dV 


2,.sin<Ju=-|^[y + F^ 


at • 


Rcbs^ 


r_L 

31 


(Hu) + ^ (Hv cos 
o? 




( 1 ) 


H = hW,X) + g. - 

The first two equations are essentially a linearized special case of the Havisr- 

Stokes equations and the third is a continuity mass conservation equation. 

Here 4> denotes latitude, X longitude, u the eastward component of velocity, 

V the northward component of velocity, f, and the components of the tidal 

X 9 

force, g the gravity acceleration, E. the mean earth radius, and m the rotational 
velocity of the Earth. The ocean depth h ( 9 , X) denotes that distance from 
the water surface to the ocean bottom when averaged over a very long time period 

T - ■ 

h-Co,!) ^ limit Y 
T-^ 0 


and ^;^C 6 ,X,t) , ^ C 9 ,X,t)_ denote the-» instantaneous elevations ■ of the ocean surface 
and ocean bottom above their time averaged values. Thus the Instantaneous upward 
Gisplacement of the ocean surface relative to the ocean bottom, the quantity 
measured by tidal gauges,' is given by (see Figure 1) 

€C^,X,t) = h(9,X,t) - h(9,>0 


(2) 


?^(9,x,t) -5j^(9,x,t) 


These equations are referred to as *^long”Wave*^ equations. In their 
derivation the ocean is regarded as a homogeneous and incompressible fluid, 
and nonlinear acceleration terras. are neglected as second order effects. The 
assumption of hydrostatic equilibrium is invoked in deriving the continuity 
equation, in that it makes the horizontal velocity components independent of 
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the radial direction and allows integration of the escpression for the balance 
of mass in a vertical water column froni the ocean floor to the ocean surface. 

The vertical velocity W need not be considered escplicitly^ as the Coriolis term 
,(i3COS(|)W is neglected in the first equation for the eastward velocity com- 
ponent. 

The displacements of the ocean and solid earth in response to T, the direct 
tidal forces due to the moon and sun, results in lxl altered mass distribution 
due to the formation of 'the tidal bulge inside the Earth and hence a change 
in the Earth’s gravitational field. Let and denote the gravitational 
disturbances of the bodily tide and ocean tide, respectively. Moreover, 
denote the forces of friction and lateral turbulent viscosity by and E^. Then 

F = T + E^ + E + F . (3) 

b s f V 

The direct tidal force, T, can be expressed as the gradient of the tidal 
potential E, 

T vr . (4) 

The gravitational disturbance due to the bodily tide results from two sources; 
a primary response to the non-loading astronomical potential and a secondary 
response to the crustal loading by the time varying ocean tide. The Green’s 
functions for the augmented potential and elevation caused by the loading 
of the crust at an angular distance Y (spherical earth) per unit of loading 
mass are given hy Farrell [10] as 

.00 

5^(T) ^ ^ -Lj k* P (cosy) 

e n=^0 ^ ^ 

■ ' ■ ■■ ■ ■ ■ e " 

where M is the mass of the earth and Ic’ , h’ are load Love numbers. The 
e n n 

primary solid earth tide, due to the, fact that the free periods of the 
elastic earth vibration are less than an hour, is accurately represented as a 
static response. Thus the total body elevation and gravitational disturbance 
axe"' 
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^ = 


^ r 

g 


R U* (y) dS2^ 


( 6 ) 


F, = k^T + W, 
b L b 


(7) 


where 


= I I R C P^Uy) d 


( 8 ) 


Here p is the mean density of sea-water, the domain of integration is over the 
global ocean surface and y is the angle between the integration variable 
and the position at which the integrals are being evaluated. 


I ; 

Li 


i i. 


cosy — sin4)sin^’ + cos0cos<t)’ cos(X-X^). 


(9) 


The first terms (static) in the equations for and are expressed in 
terms of the Love numbers and h^ with values of approximately .3 and .6, 
respectively. These parameters express the elastic yielding of the Earth 
to a non-loading potential* The gravitational disturbance of the ocean tide is 


F = VV 
s s 


(10) 


where 


V 


■Ij 


lx - xM 


( 11 ) 


1;?: 


and the domain of integration is over the ocean surface. Here G is the 
gravitational constant and x, x’ represent the position on the earth’s 
surface at which F^ is evaluated and the position of the element of inte- 
gration, respectively. On the surface of the spherical earth. 


|x - x"| - J 2R^ (1-cosy) 


j ‘ 


6 


so that 


V ^ 



GpR^^ 


2Ksin Cy/2) 


d 


£J 



R 




dS2' 


( 12 ) 


M 2sin(Y/2) . 

a 


In modeling the friction and lateral Triscostty terms we follow Zahel [25] 
and set 




/ u 


H 


:) 


R = a A 
V hv 


(;) 


(13) 


where C and C, are the coefficients of friction and hori.zontal eddy viscositv. 

r hv ’ 

respectively, and A is the horizontal Laplacian 


A “ 


2 2 
1 3 4 1 i 


2 2 2 2 2 * 
R cos f 31 R 3^ 


The Laplace tidal equations thus hecomes 


- 2wsin^v 


sr 


li „ 




Rcosi^ 31 Rcos(i> 31 r H hv Rcos^ 3X 


3v 

at 


•h 2tosin9u “ 


_ 3F ^ li ^ r 


R. 


a^f. R 3^1 . . r H 


, ^ * ,1 ^V' 

^ (14) 


where 


M , 1 [ 3 (Hu) , 3 (ttocos j>) 1 ^ . 

&t: ' 'Rcpsi!)' L ;v31 :V. . . -* 


H = h(^,>.) + 5 


viv. 




and is the potential due to oeean loading and self -gravitational effects. 


r^(<i>,A,t) -J jR'^?C4»SA*,t)p y ' 2sin( ' 72 )~^ ^*Ct) - S (15) 

Si ® 


COSY - siiWftsin^ ’ + cos^cos^ ^cos (X-X ’) 


The LTE thus becomes 5ntegro^differential in nature (IDLTE) i^hen ocean 
loading and self -gravitation effects are included, as pointed out by Hendershott. 
0.5]. . Neglecting the perturbing potential TV reduces the systenx to that of 
report A, The present study solves the IDLTE system by successive approximations 
using the time stepping method of Hansen and finite differences, imposing 
the boundary conditions of impermeable coasts, i.e . 7 the perpendicular component s- 
of tidal velocity are required to vanish at the coasts. ■ The procedure is to 
first generate a zeroth order solution of Equation (14) neglecting T ^ . Then' 

F,V is calculated from Equation (15) using the zeroth order splution for the 
tdLde A first order solution of Equation (14) is then generated using the 
potential calculated from? the zeroth order solution, . etCv. ^This successive^- ; 
iteration technique is continued until convergence is achieved.- 


Kie periodic tidal solution is represented /in the form 


where contours of equal amplitude A(^^ ,X*) are denoted . corange lines and lines 
of constant phase W are denoted co tidal lines> The expression for T^ then becomes 


F“(p, X^t) = P(<f>j,X) coscft - Q(6,X) sinot 


where 



’)cos[?(|',A')3 

1 

tn f V / 

' 

(y)-gh^(Y) 


a 

1 f 



sinijj^d^MX^ 


3,0 TIDAL gOTEHTIAL E}CPANSION 

Tlie. complex motion of the moon and stin results in a tidal disturbing force 
which displays periods ranging from thousands of years to hours. Doodson [7] 
and recently Cartwright [ 4] have expanded . the gravitational potential of the moon 
and sun at a point on the Earth surface into the series of tidal “components" 
or "constituents" 




P/Cn^jn^j ,.,ng)^p iCn^T-fn^sh, . 


( 20 ) 




2’" 



where the sma is over all positive and negative values of the integers n^. The 
arguments of the exponential are the six as trdhomlcal frequencies presented in 


Table I. The number 


(nj^jn^+S ^n^d-S jn^i*5 sngd-5) 


( 21 ) 


is knotTO as the Boodson number and uniquely identifies individual components. 


The functions $ are grouped into three species, for n^ = 0, 1, and 2. The 
functional forms are 

SPECIES I; 

3(0, , . , ,ng) = gGKP^ (sin^) 


SPECIES II: 


1 4^-y 1 

BCl.....ng) = gCKP2(sinij>)e™ 


( 22 ) 


SPECIES III: 


BC2,...,Tig) “ gGKp|(siii(!))e’'^^^ 


The CCn^sn^j • 


H m- 

.Hg) are constants, the P^ are Iiegendre polynomials and 


K - “ C- 


M 


moon. V 


2 


-) 


“ 53.7 cm 


R" 


(23) 


moon 


10 


where R denotes the mean lunar distance, 
moon - 


The fictitious displacement 


Ut) = 


3 (n- , , . . ,n^) 


g 


( 24 ) 


is denoted the "equilibrium tide” and represents a sea surface of constant 

gravitational potential. Ihe tidal constituent potential is then T = gC. 

% 

From the definition of it is clear that Species III tides are semidiurnal, 
Species II diurnal, and Species I long period tides, Doodson computed the con- 
stituents with amplitude greater than 10 ^ of the amplitude of the largest tide, 
some 400 terms. The most prominent terms are presented in Table II. 


For the present study, long period tides are not considered and for input 
purposes to the computer program, the tidal consituents are represented in the 
form 


TYEE III: 


TYPE XI: 


- CR 2. i(ot+2X) 
I - — cos ^ e ' 

■? = "CK Sin2ij) 


(25) 


where the constants C and a are specified by the particular constituent of 
interest, and 


cr ~ n^T d- n^s + “1 n^p + n^N + 


(26) 


where the rate of change is per mean solar hour. Note the factor of 2 difference 
is the semidiurnal and diurnal representation of Equation (25) . 
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Rate 

Dafitiition Period (radians per Preqnency 

mean solar hour) (deg/ solar hour) 


T 

Greenwich' celestial 

1 lunar day 

.2529359 

14.492 


loagitude 




s. 

Iiunar ; orbital 

1 tropical. 

.0095819 

,549016 


longitude 

month 



h 

Sola?: orbital 

1 tropical^ 

. .007163 

.041068 


longitude 

year 



P 

Lunar perigee 

8-85 Julian 

. 00008168 

.004642 


longitude 

years 



Hi 

Lunar node 

18.61 Julian 

.0000377 

-.002206 


longitude 

years 




Solar perigee 

20 j 900 Julian 

• 344 X 10 ^ 

.00000196 


longitude 

years 




Table I 
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Desxgaation 

Dart-rin Doodsou Kuirfcer 


Period 

(Hours) 


Amplitude Factor 




Frequency 
(radians/ day) 


Type I 



(0,7, 5,5, 5, 5) 

.156 

.459969 

327.84 

Type II 

0, 

X 

(1,4,5,5,5,5) 

.377 

5.840445 

25.82 


‘(1.6, 3, 5, 5,5) 

.176 

6.265983 

24.07 


(1,6, 5,5, 5, 5) 

.531 

6.300388 

23.93 

Type III 

"^2 

(2,4,5, 6,5,5) 

.174 

11.910864 

12.66 

“2 

(2,5,5,5,5,5) 

.908 

12.140833 

12.42 

"2 

(2, 7, 3,5, 5,5) 

.423 

12.56637061 

12.00 

’^2 

(2, 7, 5,5, 5, 5) 

.115 

12.600776 

11.97 


i I 


Tatle II 


- j. 

■; i' 
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4.0 GOMPirXATIOl^AL TECHNIQUES 
4*1 Hansen Grid 

The numerical technique employed in the present study for the Laplace Tidal 
Equations (14) follows Hansen [12] and Zahel [25] in applying a finite difference 
algorithm in temporal and spatial variables and integrating forvjard in time 
(timestspping) , We employ the global grid system proposed by Hansen (also known 
as the Richardson grid) whereby a staggered grid with velocity components and 
tidal elevations are computed at adjacent, but not coincident, grid points. The 
coastal boundary conditions are represented by horizontal and vertical lines 
imncsed on the grid and pass only through those points denoting velocity com- 
ponents. .Vertical boundary lines pass through u- grid., points (eastward velocity) 
and horizontal boundary lines pass through v grid points (northward velocity) so 
that the requirement of the vanishing .of perpendicular components of velocity 
forces, all velocity components falling on the boundary line' to be identically - • 
zero for a1 1 time. Moreover, it is to be noted that with this procedure there 
is never a corner point involved in the boundary conditions, Thiq^ is of consid- 
erable numerical importance. The Hansen grid is depicted in Figure 2. 


4.2 Finite Difference Equations 

The finite difference formulation for the Laplace Tidal Equations (14), using 
central differences for spatial variables and forward differences for the time 
variable are, following Zahel [251 

u(N, M, 

\ C At Ju(N,M,t ^[v(N,M-l,t )+v(H,M,t)+v(Nbl,K-l,t )+v(N+l,M,t )]^ 

/ _ r n iO n n n 


= 41 - 


( 


j[h (N,M) -fh (M,M-X) -i-^(N ,M, (N,H-1 , t^) ] 


u(N,M, t^) 4* 2mA t sin[4>(M)] 


. j [v(N,>J-l,t^)-t-v(H,M,t^)4'v(N+l,M-l,t^)+v(M4-l,M,t^)3 


-P C, At 
hv 


(u(K,Ji-X.t )+u(N,tK-l,t^)-2u(K,M,t ) 
) ■ n n XL 

\ [Rcos[4-(tI>]AX]^ 







Rcos[iJ^(K) ] 


+ Fjfit 


n 

AX 


|fOJ.M)cos[at^] - |fi(H,M)siii[at^] 


\il 




j CJ>t ■i[u(K-l,M,t^)+uCN-l,M!-l,t^)+u(N.lK-l,t^)+u(N,M,t^) 

) l-thCNM-O+hCK-l.iO+gCH.M.t )+g®-l,M,t ) ] 

I Z XI TX 


, v(NjM,t^) - 2tjiAt sin[A(N)- -^] 


, )4*uCR-l,yfrl,t )+uCN,K-t-l,t )+u(H,M,t )] 

<r n n n n 


Li 




( vCN,M-l,t^)4-v(N,I:-H-l,t^)-2vCN,M,t^) 


n 


n 


[Kcos[^(N)- -^3 AX 3^ 


(27) 


v(N+l,H, t^)+v(N-'l,H,t^)-2v(N,M, t^) ) 

[RA?3^ I 

(C(N,H,t )) At (ap 3o’ 

- ^ ^ ^ 3 " |§(N,M)sin[crt ]: 


K 


AO 


R 30 


n" 30 


y 


IJ 


ii 


-'■ V‘ 


■where 


K fvsvi', = ?(H,M) - P(N.M-1) 
3A AA 


1! 


M /„ Ml = Q(H,M) - Q(H.M-l) 
3A AA 


3P _ PCH.M) - P(H-l.M) 

A^ 


( 28 ) 


la f„ „■> - Q(M.M) - Q(M-1,M) 

3^ 


where 


A Rcosiji 3A 


F. = 




aacl 


SPECIES III: 


SPECIES II: 


/ =,- ^ cos I^®)] sia [2(X(M) - + ot^3 


- -^)] cos [;2X(M) +ort^3 


( \ = - siu E>^(H)] sin E(XCil) - -y) + 


-4 “ [2«(N) - ] cos [XCM) + ot_3 


rr 


(29) 


3P 30 3P 30 

Ihe quantities *^5 * 9 ^ arise from the perturbing potential and 

are held constant in time over an iteration (see Section 2.D). Por computing efficiency 
these derivatives are precalculated and held in storage for the integration 
algorithm. The continuity equation becomes 


'f 


E.cos[A(N)] I AA 

^ _ -- -j^—— . ^ 

|[h (K,M)+h (H-l®)+5 (N,M, t^)+?(N-l,Mj t^) 1 V(N,lt, t^j.p 

Ac^. 


AX/ COS (N)“ “^] 


|•[h(N-Hl,>^)'hh(N,MHeCN-hl,M,tJ4-?(N,H,t^)]v 


A^/cosf iJi(N-Kt)— ^ J 
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! 


( 30 ) 










It is to be noted that the tidal elevation ^ at time t is calculated in 
terms of the velocities at time t The finite difference formulas for the 

velocities u and v involve only nearest neighbor points, with the exception of 
the Laplacian terms for turbulent viscosity. lor grid points adjacent tc coastal 
boundaries, the evaluation of the Laplacian central difference term involves a 
grid point which would fall on land. The computer algorithm handles this by 
assuming the . normal derivative of velocity at the boundary to vanish, thereby 
eliminating the need to evaluate the point. In this way values at land grid 
points are never involved in the algorithm. (Coastal boundary points are not 
considered to be land points.) Note that these Laplacian terms are small relative 
to the other terms, having an divisor. 


4.3 Stability 

The stability of the system described by Equation (14) neglecting 

the terms has been discussed in Report A. In particular, it was shown 

that the global tides computed by the time stepping procedure become periodic 

in time and independent of the initial state (providing the initial values 

of u, were small enough that the system did not immediately blow up) 

within five tidal periods, provided the time step At is small enough to 

insure stability. The strong dissipation of the friction and lateral 

viscosity terms are the cause of this rapid convergence, quickly damping 

out the initial transients in the numerical integration and eliminating . 

spurious waves. A stability analysis for Equations less the T* terms 

has been performed by Zahel ^5] where formulas for the maximum time step 

At as a function of grid spacing, latitude limits and dissipation parameters 
max 

and are given which will guarantee stability for the tidal solution. 
Report k describes a numerical procedure for easily establishing a suitable 
time step , 

The successive approximation method proposed in Section 2-0 presents 
no difficulty in establishing a At for a given iteration to insure stability, 
as the r* terms on the right hand side of Equations 0-4) are known functions of 
position and time. The conditions under which a successive approximation 
procedure will converge to the solution of the integro-differential system 
described by Equations (l 4) perturbing potential T* will not be 

pursued in this report. 
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^ !Cha program implementation of the integration algorithm for Equations 

(14) and their boundary conditions is briefly described in Report A, The basic 

design is such that calculations are performed on a single latitude band at 

a time and required quantities are accumulated over three bands. The computer 

code features compact arrays for all variables u, v, | and h with values 

only at non-boundary ocean points on the global grid. Compact arrays also exist 
30 30 

for . The precalculation of the derivatives of P and Q increases 

the core requirements of the program, but results in a substantial savings of 
CPU time. All variables required by the integration algorithm reside in core 
vhiclt helps keep I/O and CPU time requirements down. The time needed to integrate 
Equations (L^ on a global 3 ° x 3° grid with a 3 minute step size for 12 hours 
is approximately 3 CPU minutes on an IBM 360/91. Thus, to integrate over 6 
periods of a semi-diurnal tide requires approximately 20 CPU minutes . The CPU 
time increases substantially for finer global grid meshes (and consequently 
smaller time steps). 

The program is designed with an overlay structure to minimize computer core 
storage requirements. A solution with a global grid resolution of S'* x 3 ° 
can be acconnnodated in under 50 OK bytes of core* A 2 ° x 2° global grid requires 
approximately 60 OK bytes of core, while a global 1*^ x 1“ grid would need 
slightly in excess of lOOOK bytes. 

The tidal solution 4(^,1; t) obtained by the software system in time-stepping 
from time t=0 to t=T is specified on the global grid at time t = T. It 
is convenient (and conventional) to represent the periodic solution in the 
f orm (Equation (16)) 


C(o,X;t) = A(^,A) cos [^(^,A) -b crt] 



where t is expressed in mean solar time* The computer sof tt 7 are establishes 
this representation by integrating the solution over a quarter of a tidal period 
from.t “ T to t ~ T so that 




J + g(^),X;T-hr/2ar)^ 


tan 


”1 


-g(ji.X;T+^/2a 


- OT 


(31) 


The software system employs a graphics package to display solutions in 
terms of corange lines (contours of equal amplitude A(<{> 5 A)) and cotidal lines 
(lines of equal phase I'or purposes of display, we employ the representation 


^((|»,A;t) ^ A(^,A) cos[cJt*-S(^,A)] 


(32) 


where 


S(<J.,A) - - W(^,X) . (33) 

The cotidal values appearing on the tidal charJtS presented in this report are 
expressed in hours 0 to 24, instead of angular measurement. Hour values are 
obtained by dividing the phase function S expressed in degrees by 15 degrees per 
solar hour. Note that both diurnal and semi-diurnal tides have cotidal values 
expressed in hours from 0 to 24, These relationships are given explicitly as 
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4.5 Evaluation of Greenes Functions for Ocean Loading 

The Greenes functions for the augmented potential and elevation caused 
by the loading of the crust given by Equations (5) are slowly converging 
series, as discussed by Earrell M - Eollowing Earrell, and using the 
identities 



. E^(cost). 


2sin(y/2) 


P^Ccosy) _ 
n-0 n 


£fi|2sin(y/2) j(l+sinCy/2)J | 


( 34 ) 


, we rewrite'. Equations (5) in the form 

( 


u' (y) = ^ I ^ ^ ^ 


R 


M 


e ( 2sin(y/2) n=0 


+ S K-^co ) \ } 


( 35 ) 


(y) = 


„ ( ~ /nk^ ^ \ ^ (cosy) '1 

- ? £ V~« * s. E 

e I n=0 n=l 1 


e I n==l 


nk’ - k 
n 0 ^ 

n 


P^(cosy) - k ZVL ^2sin(y/2) [1 + sin(y/2)] 


where 



= lim 


n-^ 


h» 

n 


“ 1 nk ' 

n->™ n 


( 36 ) 


We have investigated several summation techniques for evaluating the series 
in Equation (35) > including the Shanks convergence acceleration method [23]. 
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Howeverj we found the Euler txansfonnation method discussed by Earrell to be 

the most efficient, The values of and h* used are those given by Farrell 

iL n 

for the G”B Earth Model, where 

k - 2,482 

03 

h ^ - 5.005 . 

CO • 

Moreover;, we have used the value 

ho ^ " .134 

obtained by Longman [1^ for a Gutenberg Earth Model, Farrell gives values of 

h’ , k^ for 
n n ' 

n - 1,2,3,4,5,6,8,10,18,32,56,100,180,325,550,1000,1800,3000,10000. 


Quadratic polynominals were fitted through these values to obtain k^, for 
other values of n required for the summations. 

Table III presents calculated values for h^(y), ®^(y) and the Greens function 




^ 2S14/2) 


( 37 ) 


12 

For convenience of comparison, the values of and are scaled by (Ry) x 10 , 

with y in radians. The quantities are in mks units, and an applied load of 

1 kg. is assumed. Values taken for R, p, g and M are 6371000., 1030., 9.82, 

24 ^ 

and 5,975 x 10 , respectively. 


TABLE IXi 


U’ (V) K R X Y X IQ' 

'"13.19-6V. 

- 9 v 725 : 


X K;X:Y X 10' 


0,293 

0 . 838 . 

1,674 

2.076 


1.647 
0.924 
^ 0.029 

- 1.138 

- 2.311 

- 3.439 

- 4.506 

- 5.438 

- 6.263 

- 6.660 


~ 1,663 
- 1.15 9 . 
- 0 . 025 ; 
" 1.368 

2 > 886 
: 4.494 
6^1581 
7.850 
9.580 
11.367 
31,214 
53.987 
78.671 
104.017 
128.809 
152.136 


192.523 

209.427 

224,447 

238.011 

250,661 

262.991 

275.579 

288.956 

303.668 

320.324 


Green ’ s Furic tl on 
(m^/sec\g) ^ . 
73.178 
30.777 . 

17 , 576 '' ■ . 

11.964 ■ ■ 

9 vl 44 • 

.., 7,310 
6.112 
5.339 
4.798 
4,311 
2,328 
1.716 
1.406 
1.238 
1.142 
1.089 
1.061 
1.046 
1,039 
1.035 
1.033 
1.032 
1.030 
1.028 
1.027 
1.027 
1.025 


Solutions for the principal tidal constituents S^j 1^2* ^ 2 ’ ^1’ *^1 


and have been generated using the integration algorithm described in 


Section 4.1 neglecting the terms involving T* and are displayed in Figures 3-9* 
These tide models are then solutions of the LTE as described in Report A, 
and do not exhibit an integro-dif f erential character (IDLTE) . They would, 
however, serve as the initial approximations for the successive approximation 
algorithm described in Section 2.0 for the IDTLE (see Section 4.7). Global 
grid resolutions of both 3° x 3® and 2® x 2° have been used for generating 
tide models for the principal constituents and display excellent agreement. The 


Love parameter values used in these solutions were = ,3 and h^ = ,6 while 
the dissipation parameters were = .003 and = 10^, respectively (the 
same values used hy Zahel [25, 26]). In the 2° x 2® numerical models the 
effects of the tides of polar seas north of 80° latitude are ignored 


(81° latitude for the 3° x 3° models) . This is due to the fact that increasing 
the northern latitude limit necessitates decreasing the integration step 
size to maintain stability, hence greatly increasing the computer time needed 
to arrive at the solution. However, our ntimerical experimentation has shown 
that increasing the northern latitude limit has no discernable Impact on the 
solution south of 80°, The ocean depth profile used in the solutuions is 
presented in Figure 10. Initial values of u, v and? at time t = 0 were set 
to zero. A time series plot of the tidal solutuion, ?, for a grid point 
in the central Atlantic is presented in Figure 11 . The solution grows from 
the zero initial solution to a steady state solution within several periods. 

All solutions for the principal constituents presented in this section were 
integrated over eight tidal periods (assuring a steady state solution) before the 
amplitudes and phases (Equation (31) ) were calculated. 


Prominent features evident in both the diurnal and semi-diurnal solutions 
are amphidromic systems, or geographic regions where the tidal amplitude 
vanishes. Cotidal lines radiate from the amphidromic points, where the sense of 
the amphidromic system is denoted as the direction in which the high tide 
occurs progressively later, Amphidromes are usually (not always) counter- 
clockwise in the Northern Hemisphere and clockwise in the Southern Hemisphere. 
These stationary centers can result physically as a consequence of the character- 
istic wave velocity being of the same order of magnitude of the velocity of 
the sublunar point. 
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The family of semi-diurnal waves exhibits very similar character- 

istics in the tidal solutions, as dpes the family K^, 0^, of diurnal waves. 

¥or this reason we wixl only discus'^ the and solutions of Figures 3 and 
7 • in detail* \ 

The K 2 solution of Figure 3 Was compared with the numerical., solutions of 
Hendershott [151, Peheris and Accad [21], Zahel [25] j and Bogdanov and Magarifc [3] as 
to tidal amplitude, phases and locatn^on of amphidromic systems -These solutions are con 
-veniently siimmarized hy Hendershott [14], Our solution compares very closely 
with Zahel, ^fith nearly a systematic 30% decrease in tidal amplitudes consistent 
with the factor (l-l-k^-h^) in our solution. The various models are in fairly good 
agraenient in the North and South Atlantic, with the exception that Hendershott 
does not find the South Atlantic amphidrome, and the amphl'drome at the southern 
tip of South America is farther west in ours and Zahel ^s solution relative to the 
dissipative model of Pekeris. All models are again quite similar in the Indian 
Ocean, being dominated by a large ”antiaraphidrome*' , i.e. a region of high tidal 
amplitude moving nearly synchronously. The amplitude of this feature sho’/m by 
Ker.dershott is substantially larger than the other models- Our experimentation 
shows that the region at the southern tip of Madagascar is near an amphidromic 
state. The Pacific basin shows the greatest differences between the numerical 
models. Our solution, along with Zahel, shows the appearance of an amphidromic 
system off of California, and a line of rapid phase change in a concave arc 
running from Central America to Southern Pern. Our model, in agreement mth Zahel, 
shcv.s two antiamphidromes on either side of the dominant South Central Pacific 
-amphidrome, although the eastern antiamphidrome of our model Is more pronounced 
than that of Zahel, The Western Pacific -..basin does not show the pseudo^aodal 
line extending from Japan to New Guinea present in the empirical chart of 
Dietrich [6], This region in our model is instead dominated by an amphi- 
dromic system. The full 360“ counter clockwise phase change about New Zealand 
is in agreement with the empirical charts , 


u 

u 


n 

u 


I 


n 


The tidal solution in Figure 7 was compared with that of Zahel [2Q . 

The amphidromic systems of the Atlantic and Indian Oceans correspond quite closely 
between the two solutions, with the exception that the aero co tidal contour 
runs directly between the North Atlantic and Central Atlantic amphidromes of 
our model rather than to the coasts of Spain and Dakar, respectively. As with the 
tide, the Pacific basin offers the most contrast between the numerical models. 
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Where Zahel finds a single dominate amphidrome in the Central Pacific, we find 
this region divided into two centers along an arc running from Wake Island 
to Tahiti* In addition we find an amphidromic system off the coast of Chile, 
two models compare closely in the northeast Pacific and the South Pacific. 

As with the wave, general agreement is found with the empirical chart of 
Dietrich with the exception of the western and southeastern Pacific hasin. 

In particular, the region of constant phase in the western Pacific found by 
Dietrich is dominated by an amphidromic system in our model. 

Tables IV-V compare tidal amplitudes and phases of the 
models at selected locations with the harmonic constants available from the 
Admiralty Tide Tables [1] . The phases are relative to the Greenwich meridian, 
and are the S (^ , y) phases defined hy Equation (33 • 
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Figure 4 

2® X 2® (Neglecting Loading and Self-Gravitation Effects) 

S: Dashed lines cotidal lines (phase) in hours 

A; Solid lines corange lines (amplitude) in decimeters 

C(4i,^;t) = A((|),X) cos[ at - S] 

a = .0001454452 radians/sec 







Figure 5 

2® X 2® N 2 (Neglecting Loading and Self-Gravitation Effects) 

S: Dashed lines co tidal lines (phase) in hours 

A; Solid lines corange lines (amplitude) in centimeters 

^((J),X;t) = A((|),X) cos [at - S] 

a = .00013788 radians/sec 
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Figure 7 

2° X 2“ (Neglecting Self-Gravitation and Ocean Loading Effects) 

S: Dashed lines co tidal lines (phase) in hours 

A: Solid lines corange lines (amplitude) in decimeters . 

C(^,X;t) = A(<^,X) cos [at - S] 

a = .0000729216 radians/sec 
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Figure 9 

2® X 2® (Neglecting Self-Gravitation and Ocean Loading Effects) 

S: Dashed lines cotidal lines (phase) in hours 

A: Solid lines corange lines (amplitude) in decimeters 

^((J),A;t) = A((J),X) cos [at - S] 

0 = .0000725236 radians/sec 
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Figure 10 

GLOBAL OCEAN DEPTH PROFILE 


Contour Values Expressed in Hundreds of Meters 
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TABLE IV 

Comparison of and Tide Models with Observations 
(Neglecting Loading and Self-Gravitation Effects) 



M^ OBSERVED 

MODEL 

OBSERVED 

S^ MODEL 

Azores (N. Atlantic) 

63.3° 

26.3° 

82°. 5 

48° 

(37° A4’ N. 25° 40* W) 

49 cm. 

39 

18 cm. 

18 

Bermuda (N. Atlantic) 

20° 

300° 

45° 

310° 

(32° 19’ N. 64° 50’ W) 

38cm. 

19 

8 cm. 

8 

St. Helena (S, Atlantic) 

87° 

54° 


70° 

(15° 55’ S. 5° 42’ W) 

34 cm. 

25.5 


11 

Chagos Is. (Indian Ocean) 

270° 

252° 

300° 

268° 

(7° 21’ S. 72° 27’ E) 

52 cm. 

26 

27 cm. 

13 

Cocos Is. (Indian Ocean) 

306° 

266° 

5° 

288° 

(12° 5’ S. 96° 51’ E) 

27 cm. 

22 

9 cm. 

12 

San Francisco (E. Pacific) 

215° 

195° 

210° 

200° 

(37° 48' N. 122° 27* W) 

54 cm. 

35 cm. 

12 cm. 

14 

Honolulu (Pacific) 

61° 

35 

62° 

63° 

(21° 18’ N. 157° 52’ W) 

16 cm. 

33 

5 cm. 

13 

Easter Is. ( S, Pacific) 

32° 

281° 

240° 

300° 

(27° 9’ S. 109° 27’ W) 

20 cm. 

35 

6 cm. 

16 

Apia, Samoa (S. Pacific) 

171° 

87° 

150° 

105° 

(13° 48' S. 171° 46’ W) 

38 cm. 

54 

9 cm. 

24 

Siapan (S. Pacific) 

300° 

159° 

315° 

190° 

(15° 12’ N. 145° 43’ E) 

17 cm. 

21 

5 cm. 

9 

Yaruto, Marshall Is. (S, 
Pacific) 

125° 

88° 

150° 

109° 

(5° 55’ N. 169° 39’ E) 

47 cm. 

46 

26 cm. 

23 

Burgana, Solomon Is. (S. 

Pacific) 

134° 

196° 

180° 

210° 

(9° 11’ S. 160° 13' E) 

11 cm. 

33 

8 cm. 

15 


TABLE V 


Comparison of and 0^ Tide Models with Observations 
(Neglecting Loading and Self-Gravitation Effects) 



OBSERVED 

MODEL 

0^ OBSERVED 


0^ MODEL 

Azores (N. Atlantic) 

66® 

70® 

315® 

• 

80® 

(37® 44' N. 25® 40' W) 

4 cm. 

4.6 

3 cm. 


2 

Bermuda (N. Atlantic) 

195 

230® 

210® 


220® 

(32® 19' N. 64® 50’ W) 

7 cm. 

4 

5 cm. 


3.5 

St. Helena (S. Atlantic) 


165® 

- 


137® 

(15® 55’ S. 5® 42’ W) 


3 



2 

Chagos Is. (Indian Ocean) 

84® 

80® 

60® 


62® 

(7® 21’ S. 72® 27' E) 

3 cm. 

4 

3 cm. 


3 

Cocos Is. (Indian Ocean) 

165® 

168® 

127® 


139® 

(12® 5’ S. 96® 51’ E) 

11 cm. 

11 

9 cm. 


8.5 

San Francisco (E. Pacific) 

228® 

229® 

210° 


193® 

(37® 48’ N. 122® 27’ VJ) 

37 cm. 

43 

23 cm. 


29 

Honolulu (Pacific) 

225® 

266® 

225° 


230® 

(21® 18’ N. 157® 52’ W) 

15 cm. 

19 

8 cm. 


15 

Easter Is. ( S. Pacific) 

18® 

57® 

9® 


17.5® 

(27® 9’ S. 109® 27’ W) 

8 cm. 

8 

5 cm. 


8 

Apia, Samoa (S. Pacific) 

75® 

28® 

60® 


340® 

o 

(13® 48’ S. 171® 46’ W) 

5 cm. 

7 

3 cm. 


3 

Siapan (S. Pacific) 

75° 

130® 

45® 


82® 

(15® 12' N. 145® 43' E) 

14 cm. 

11 

10 cm. 


11 

Yaruto, Marshall Is. (S. 

75® 

306® 

45® 



Pacific) 






(5® 55' A, 169® 39' E) 

9 cm. 

3 

6 cm. 



Burgana, Solomon Is. (S. 

45® 

60® 

22® 


57® 

Pacific) 






(9® 11' S. 160® 13’ E) 

20 cm. 

8 

11 cm. 
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4.7 


M 2 Tide Including Ocean Loading and Self-Gravitation 


The solutions of Section 4.6 show general agreement with empirical solutions 
with respect to the positions of amphidromic systems and tidal amplitudes. 

Notable exceptions to the agreement of the theoretical with empirical solutions 
are in the western and southeastern portions of the Pacific basin. Phase relations 
also differ in portions of the North and South Atlantic. To improve the M 2 
model of Section 4.6, the method of successive approximations outlined in 
Section 2.0 was used to solve the IDLTE given by Equation (L4) on a 3° x 3° 
global grid. The M 2 solution of Section 4.6 was used as the initial solution. 

The convergence of the procedure was slow in the southeastern Pacific 
basin, requiring 16 iterations to reach final convergence. Convergence in all 
other regions was more rapid, requiring approximately 8 iterations, (con- 
vergence was considered to be achieved when the corange and cotidal lines 
no longer changed from one iteration to another) . To assure that the iteration 
procedure had stabilized, two additional iterations (18 iterations) were 
calcualted. Figure 12 displays the self- gravitation and ocean loading 
potential F’ obtained by evaluating Equations (15) with the initial solution 
from Section 4.6 ( the 0th order approximation). The portions of F* due to 
self-gravitation only and crustal deformation only (U^ Green’s Function) are 
displayed in Figures 13 and 14 respectively . The converged M 2 solution (16th 
iteration) is presented in Figure 15 . 

Comparison shows no large differences in amplitude or position of amphidromes 
between the M 2 solution of Section 4.6 and that of Figure 15 which includes 
the perturbing effects of ocean self-gravitation and crustal loading. There are, 
however, significant changes in phase. In particular, the phases of the North 
Atlantic amphidrome are rotated clockwise approximately 25® and the phases 
from the South Atlantic amphidrome are rotated clockwise by 45° along the 
South American coast and approximately 20® along the equatorial coast of 
Africa. These shifts in phase improve agreement with observations in the 
Atlantic as demonstrated by the comparisons at Bermuda, the Azores and St. 

Helena in Table VI • The rotation of phases by approximately 20® in a counter 
clockwise manner about the amphidrome off the west coast of Australia likewise 
improves agreement with observations in the Indian Ocean, as demonstrated by 
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Figure 12 


3® K 3° M_ Potential — (ON ZEROTH ITERATION) 
2 g 


S: Dashed lines cotidal lines (phase) in hours 

A: Solid lines corange lines (amplitude) in centimeters 
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Figure 13 


3® X 3® M 2 Self-Gravitation Potential only 


(ON ZEROTH ITERATION) 


S: Dashed lines - - - phase lines in hours 

A: Solid lines amplitude lines in centimeters 





Figure 14 


3® X 3® Crustal Response to Ocean Tidal Loading (ON ZEROTH ITERATION) 

S: Dashed lines phase lines in hours 

A: Solid lines amplitude lines in centimeters 
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Figure 15 

3® X 3° M 2 (Including Loading and Self-Gravitation Effects) 

S: Dashed lines co tidal lines (phase) in hours 

A: Solid lines corange lines (amplitude) in decimeters 

C(4>,^;t) = A((J),X) cos [at - S] 

a = .00014052 radians/sec 






TABLli Vi 


Azores (N. Atlantic) 

(37° 44* N. 25° 40' W) 

Bermuda (N. Atlantic) 

(32° 19* N. 64° 50' W) 

St, Helena (S, Atlantic) 

(15° 55' S. 5° 42’ W) 

Chagos Is, (Indian Ocean) 

(7° 21' S, 72° 27’ E) 

Cocos Is, (Indian Ocean) 

(12° 5' S, 96° 51' E) 

San Francisco (E. Pacific) 

(37° 48' N. 122° 27' W) 

Honolulu (Pacific) 

(21° 18' N, 157° 52' W) 

Easter Is, (S, Pacific) 

(27° 9' S. 109° 27' W) 

Apia, Samoa (S, Pacific) 

(13° 48' S. 171° 46' W) 

Siapan (S. Pacific) 

(15° 12' N, 145° 43' E) 

Yaruto, Marshall Is, (S. Pacific) 
(5° 55' N. 169° 39' E) 


OCEAN 

TIDE PAHSE AND Al-EPLITUDE 

OBSERVED 

THEORETICAL 


(Without Ocean Loading 
and Self-Gravitation) 

63.3° 

26.3° 

49 cm. 

39 

20° 

300° 

38 cm. 

19 

87° 

54° 

34 cm. 

25.5 

270° 

252° 

52 cm. 

26 

306° 

266° 

27 cm. 

22 

215° 

195° 

54 cm. 

35 cm. 

61° 

35 

16 cm. 

33 

32° 

281° 

20 cm. 

35 

171° 

87° 

38 cm. 

54 

300° 

159° 

17 cm. 

21 

125° 

88° 

47 cm. 

46 

134° 

196° 

11 cm. 

33 


THEORETICAL 

(With Ocean Loading 
and Self-Gravitation) 

50° 

41 

308° 

17 

73° 

26,2 

269° 

31 

285° 

29 

210 ° 

45 cm. 

66 ° 

31 

303° 

34 

112 ° 

55 

200 ° 

25 

114° 

53 

210 ° 

34 


Burgana, Solomon Is. (S. Pacific) 
(9° 11' S. 160° 13' E) 


comparisons at Chagos Island and Cocos Island. In the Pacific northeast, 
the phase lines radiating from the amphidrome off the California coast are 
rotated clockwise by approximately 15® improving the model along the North 
American western seaboard. The only appreciable shift of an amphidrome involves 
the system off the coast of Equador, which moves southward by approximately 
15° of latitude as well as being rotated clockwise by 45°. The amphidrome which 
dominates the northwest Pacific basin experiences a 45° counter clockwise 
shift of the phase lines along the Asian coast making phases along Japan 
and Okinawa more realistic. Agreement of the M 2 solution including ocean 
loading and self- gravitation effects with island observations is generally 
within 30° of phase and 25 cm. is amplitude. Exceptions are along the Aleutions, 
in the Sargasso Sea region and the southeastern Pacific basin. (see Table VI). 

Solutions for and S 2 including the effects of ocean loading and self- 
gravitation are currently being calculated and will be presented in a later 
publication. 

4.8 Total Tide Comparison with Observations 

The complexity of the total tide becomes evident when the seven constituent 
models are combined to yield an approximation to the total ocean tide. The 
theoretical global ocean tide models are not expected to be realistic along 
highly irregular coast lines or in shallow water along continental coasts due 

to the assumptions imposed in the boundary conditions. Deep ocean islands, 
such as the principal ports of Honolulu, Bermuda and the Azores of the Admiralty 
Tide Tables provide a reasonable means of comparing the ocean tide models 
with observations. It must be remarked, however, that groups of islands, 
to small to be included in the global model, can strongly influence the tides in 
their vicinity. Figures 16 - 19 display the tide composed of the S 2 > N 2 > 

K 2 » K^, Pj^ models of Section 4.6 and the M 2 model of Section 4.7 for 
the Azores, Bermuda, San Francisco and Honolulu. Times and values of high 

and low water from the Admiralty Tide Tables are indicated on the plots. 

The comparison at the Azores, which is dominated by the M 2 Tide, is very 
good, while the predicted tide at Bermuda displays a reasonable phase but 
small amplitude (see also Table VI) . The tides at San Francisco and Honolulu 
show a greater diurnal influence. The phases and amplitudes of the predicted 
tide at San Francisco are reasonable, while the predicted tide at Honolulu shows 
poorer agreement. Causes for possible disagreement between observations and 
the theoretical models were discussed in Section 1.0, 
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Meters 


Figure 16 


I 

i 

Tide Model Prediction for Azores, Jan. 1, 2 1977 
(Circled Points are Observed High and Low* Waters) 
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I Figure 17 

‘ Tide Model Prediction for Bermuda, Jan. 1, 2 1977 

i (Circled Points are Observed High and Low Waters) 








5*0 Goniputer Software System 

The integro- differential IDLTE Software System is composed of three 
basic modules: 

I. IDLTE Integration Package. This program itegrates the system 
described by Equation (14) of Section 2.0, It requires as 
input the prescribed continental boundaries, the global 

ocean depth profile defined on a grid, an initial solution, and 
the potential To calcualte the zeroth order solution both 

the initial solution and F* are set to zero. 

II. Surface Integration Package. This program evaluates the potential 
F’ of Equation (15) of Section 2.0. It requires as input the 
prescribed continental boundaries, a {Specified global tide 
^((}),A) in the form of gridded amplitude and phase values and a 
table of values for the Green’s functions®’ and U’ of Equation 

(5) as a function of angle. 

III. Ocean Loading Green’s Function Package, This program evaluates 
a table of the Green’s functions and U’ of Equation (35) from 
an input set of load Love numbers k’ and h’ . 

The three modules are written in FORTEAN IV, and Modules I and II 
are designed with a modified graphics package (WOLF PLOT PACKAGE) for solution 
display. The I/O design of the system uses multiple magnetic tape drives. 
Table VIII defines the I/O units for Modules I and XI, 

5.1 Program Input 

User control of Modules X and II is set up by card input for the 
principal variables through a HAMELIST specification. The variable names, 
their function and their default values are presented in Tables IX and X . The 
ocean depth values required for Module I are read from an external tape, with 
values defined on a world grid, (depth values on ocean grid points and 
0 on land grid points) . 

Following the NAMELIST INPUT for both Modiiles I and II, the continental 
boundary data is input. The ocean boundaries are instituted as horizontal and 
vertical segments on the grid with horizontal boundaries along v velocity 


points and vertical boundaries along u velocity points. The coastal boundary 
configuration is input by cards in a special compact format. For each 
latitude on the world integration grid there is a data card of the format 
14X5 with fields. 


lAT NO IDATAC12) 

where LAT is the latitude in degrees, NO is an integer 0 or 1 depending on 
whether the point at zero longitude at the particular latitude is on land or 
ocean, respectively, and the array IDATA is the value of the longitude, in 
degrees, at each grid point for that latitude where a land-water boundary 
occurs. With this input technique any arbitrary coastal boundary configuration 
-may be implemented. The boundary configuration presently employed is a 
realistic 1° x 1° world ocean boundary. The user specified grid size extracts 
the proper ocean boundary from this 1° x 1° base. 

For Module XI, the ocean boundary data is followed by a table of values 

for the Greenes Functions II* (y) and^* (y) for y=0 » 2°, ISO'^ 

-12 

(scaled by (Ry) x 10 as discussed in Section 4,5) . The input format is 
8F10,3. 


Module III is a stand-alone program designed to calculate the Green *s 

Functions and IT* as described in Section 4.5,^^® program need only be 

run when a new set of load Love number h* and k* are employed to calculate 

n n 

new functions and U* . The data for the program consists entirely of values 
of h , k and h*, k* for a set of values of n. This data is input in data 
statements within the main program of Module III, 
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•Table VII ’ 
MODULE I SUBRODTIMES 


Subroutine Katue Function 

SPIASH The ^executive program which reads Mamelist Input and calls the 

working subroutines of the system. This subroutine contains 
a BLOCK DATA subroutine. 


I-IAINA The nodule which reads the input coastal boundary data and sets 

, up the world function (0 on land points and 1 on ocean points) 

corresponding to the selected integration grid size. Prints 
the world function to display boundaries and defines the trans- 
formations from the global grid to the compacted ocean arrays. 

XIDFLT The platting module which generates corange contours, cotidal 

lines and time series plots for selected geographic points. 

FUNS The function subroutine which evaluates the derivatives of the 

tidal component for use in the driving terms of the Laplace 
Tidal Equations. 

ALGR The integration algorithm for the Laplace Tidal Equations, 


DEPTH 


SVmiTR 


The nodule x^^hich transforms the input ocean depths defined on 
a world grid to the compacted ocean grid. 

INPUT /OUIPUr routine for the world grid ocean depth values. 


ZARY 


The nodule which transforms the compacted ocean arrays onto 
world grids for display and plotting purposes. 


Table VIII 

TAPE I/O UNITS FOR MODULES I & II 


Module 

I/O Unit 

Function 

Module I 

30 

INPUT - World ocean depth profile 


40 

INPUT - Initial tide solution 


’ 50 

OUTPUT Tide solution co tidal and corange 



values 


60 

OUTPUT - Time series of tide heights over 



the span of integration for three 


70 

specified geographic location, 
ax’ a<i)’ 3X’ 3<f> 


80 

OUTPUT - Tide solution ^ at the final integration 



time. 


Module II 

30 

INPUT - Tide solution co tidal and corange 
values (from unit 50 of Module I) 


40 

OUTPUT - Tide perturbation potential T ^ 

in CO tidal and corange form. 

3P aP 3Q IQ 

OUTPUT - 


50 


52 


Table ix 


MODULE I namelist VARIABLES 


Variable 

N^nie 


Function 


Default Value 


DDELT 

TTis 

BR 

SR 

AH 

K 


K2 

CT 

SIO 

G 

W 

m 


Integration time step in minutes . 

Integration start time in minutes when the restart 
option is used. 

Radius of the earth in meters. 


Friction coefficient, 0 . 

’ r ■ 


( 6 .) 

( 0 .) 

(6371 X 10^) 
(.003) 


Coefficient of lateral turbulent viscosity, in. (1 x 10 ) 

meters /sec. 


Dimensionless ’’general lunar coefficient” 


(.842 X 10^^) 


, M 

3 moon 




2 M 


earth 


R 


moon 


where M denotes mass, R^ the earth’s radius and 

R the mean lunar distance, 
moon 

Love number which represents the coefficient of 
additional potential due to deformation by a 
non-loading potential. 


( 0 .) 


(.1405 X 10 


Amplitude factor for tidal constituent. (C of Eq (25)) (!■) 

Tidal constituent frequency in rad/sec. 

2 

Gravitational acceleration in meters/sec . (9.81) 

Rotational rate of the earth in rad/sec. (.72722 F.-4) 

Grid size in degrees. (6°) 


IFLGR Restart Flag (0) 

0 - Start integration from t-0 from zero initial 

solution, 

1 - Read tidal solution at time TTN from I/O Unit 

40 and integrate over 1/4 tidal period; Form 
corange and co tidal solutions. 

2 - Read tidal solution at time TTN from I/O Unit 

40 and integrate over RITR tidal periods. 
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Table IX (Coa t.’) 
MODULE 1 MAJyfELIST VARIA:BLES 


Variable 


U ame 

Function 

Default Value 

XPERT 

Perturbation Flag 

cu; 


0 - No perturbation potential included in 

equations of mjtion. 



1 - Read ||s ||, |S, |S from I/O Unit 70, and 



include F * potential in equations of motion* 


WOO 

Lower latitude limit in degrees for global inte- 
gration. NOTE: The difference between -90'* and 
NOO must be an integer times the grid size ND. 

(-78®) 

NEE 

Upper latitude limit in degrees for global inte- 
gration. NOTE: The difference between 90'* and 

NEE must be an integer times the grid size ND. 

(84") 

NITR 

The number of tidal constituent periods over which 
the integration is to proceed. 

(1.) 

Ni 

Grid point number (in compacted form) to be sampled 
for a time series analysis. 

(1.) 

l\2 

Grid point number (in compacted form) to be sampled 
fcr a time series analysis. 

(2.) 

N3 

Grid point number (in compacted form) to be sampled 
for a time series analysis. 

(3.) 

ITTPE 

Species of tidal component being integrated. 

2 - diurnal 

3 - semi-diurnal 

(3.) 


H2 


Love number which represents the coefficient of 
uplift response of the solid earth to a non-loading 
potential. 


( 0 .) 


Table X 

MODULE II NAMELIST VARIABLES 


Variable 

Name 


BR 


KOO 


NEE 


GFCNG 


GFHP 


GFKP 


GFSG 


Function 


Default Value 


3 

Radius of the earth in meters (should be the (6371x10 ) 

same value used to generate the tidal solution 
from Module I. 

'Lower latitude limit in degrees for global integration. (-78'’) 
(Must be the same value used to generate the 
tidal solution from Module I. 

Upper latitude limit in degrees for global 
integration. (Must be the same value used to 
generate the tidal solution from Module I. 

Greens Function Flag 

0. - Use default Greens funtion 

(See Table III) 

1. “ Use input values of and U^ functions 

to calculate Greens function. 

U' Input Flag (Valid ONLY IF GFING =1.) (1.) 

0. ” Set to zero in the Greens Function 

1, - Input U* as a function of angle (0°~180'’) 

in increments of 1'’ [F0RMAT 8F10.3] in units 

U' (G)*A*0*1o12, 

§* Input Flag (Valid ONLY IF GFCNG =1.) (1.) 

0. - Set 0* to zero in the Greens Function 

1. “ Input o’ as a function of angle (0° - 180°) 

in increments of 1° [F0JRMAT SF10.3] in units 
a’ (0) * A * 0 * lO^^. 

Self-Gravitation Term Flag (1.) 

0, - Set self -gravitation term to zero in the 

Greens Function. 

1. - Include self-gravitation term in the 

Greens Function. 


C84°) 

( 0 .) 


5. 2 Program Output 


The data output from the integration algorithm over the ocean grid network 
is voluminous and difficult to analyze in taw form. This necessitates plotted 
displays for visual analysis of the solution. In particulars the principal outputs 
are the cotidal and corange solution contours and the tidal time series plots for 
selected geographic points which are output to a plot tape for interface with an 
external plotting device. In addition, the integrator time t and the variables 
Uj V, C and h+| are output to tape in their compacted ocean arrays at the end of 
each integration period for program restart purposes. Additionally, the tidal 
time series is output on tape at the termination of the integration. 

Printed plots are also produced for every plot output to the external plotting 
device. In additions ' the program prints the total world function and the world 
grid of the ocean depth values for user inspection prior to calling the integra- 
tion algorithm. At the end of each integration period the variables u, v and ^ 
are printed in their compacted ocean arrays, and then are printed in a full world 
grid array display. 


6.0 


NEW TECHNOLOGY 


The effort under this contract consisted of the development and 
progranming of techniques to numerically calculate global semidiurnal and 
diurnal ocean tide models by integration of the Laplace Tidal Equations in 
th.e time domain incorporating friction, lateral viscosity terms and 
crustal loading effects. Global K^s 0^ and solutions 

have been obtained (neglecting self-gravitation and crustal loading) and 
compared with solutions existing in the literature. An solution including 
the effects of self-gravitation and crustal loading has been obtained by successive 
approximations and compared with observations. 

Frequent reviews and a final survey for new technology were performed. 

It is believed that the mathematical and programming techniques and 
algorithms developed do not represent "reportable items, " or patentable 
items, within the meaning of the New Technology Clause. Our reviews and 
final survey found no other items trhich could be considered reportable 
items under the Kew Technology Clause. 
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